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Popular Summary 

New state of the art methodology is described to analyze AIRS/AMSU/HSB data in the 
nresence of multiple cloud formations. The methodology forms the basis for the A 
Ic“ea7algcri,h m which will be used to analyze AIRS/AMSU/HSB data on EOS 
Aqua. Results are shown for AIRS Science Team simulation studies with multiple cloud 
formations. These simulation studies imply that clear column radiances can be 
reconstructed under partial cloud cover with an accuracy comparable to single spot 
channel noise in the temperature and water vapor sounding regions, temperature 
soundings can be produced under partial cloud cover with RMS errors on the order of or 
better than, 1°K in 1 km thick layers from the surface to 700 mb, 1 Ion layers from 
mb to 300 mb 3 km layers from 300 mb to 30 mb, and 5 km layers from 30 mb to 
and moisture profiles can be obtained with an accuracy better than 20% absolute errors in 
1 km layers from the surface to nearly 200 mb. 
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ABSTRACT 

New state of the art methodology is described to 
analyze AIRS/AMSU/HSB data in the presence of 
multiple cloud formations. The methodology forms 
the basis for the AIRS Science Team algorithm which 
will be used to analyze AIRS/AMSU/HSB data on 
EOS Aqua. The cloud clearing methodology requires 
no knowledge of the spectral properties of the clouds. 
The basic retrieval methodology is general and 
extracts the maximum information from the radiances, 
consistent with the channel noise covariance matrix. 
The retrieval methodology minimizes the dependence 
of the solution on the first guess field and the first 
guess error characteristics. Results are shown for 
AIRS Science Team simulation studies with multiple 
cloud formations. These simulation studies imply that 
clear column radiances can be reconstructed under 
partial cloud cover with an accuracy comparable to 
single spot channel noise in the temperature and water 
vapor sounding regions, temperature soundings can be 
produced under partial cloud cover with RMS errors 
on the order of, or better than, 1°K in 1 km thick 
layers from the surface to 700 mb, 1 km layers from 
700 mb to 300 mb, 3 km layers from 300 mb to 30 mb, 
and 5 km layers from 30 mb to 1 mb, and moisture 
profiles can be obtained with an accuracy better than 
20 % absolute errors in 1 km layers from the surface to 
nearly 200 mb. 

1. INTRODUCTION 

AIRS (Atmospheric Infrared Sounder) is a high spectral 
resolution (v/Av ~ 1200) infrared sounder, with 2378 

channels covering the spectral domain 650 cnT* - 2675 
cm'l, which will fly on the EOS Aqua platform in 
2002, accompanied by the AMSU A (Advanced 


Microwave Sounding Unit A) and HSB (Humidity 
Sounder for Brazil, which is similar to AMSU B). The 
AIRS footprint is 13 km at nadir, as is the HSB 
footprint, with a 3x3 array of AIRS and HSB footprints 
falling into a single AMSU A footprint. Characteristics 
of the AIRS instrument are given in Aumann et al , 
2002 . 

Susskind et al ., 1998 described the first version of the 
methodology used by the AIRS Science team to analyze 
AIRS/AMSU/HSB data in the presence of clouds to 
determine surface skin temperature, surface spectral 
emissivity and bi-directional reflectance, atmospheric 
temperature-moisture-ozone profile, and the heights and 
amounts of different layers of clouds in die fields of 
view. Two important characteristics of the basic 
retrieval methodology are that no assumptions are 
needed about the spectral properties of die clouds and 
no assumptions are needed about the intrinsic accuracy 
of the first guess field used to start die iterative process. 
This paper describes further dieoretical improvements 
in the retrieval and cloud clearing mediodology 
incorporated in the current version of the AIRS Science 
team algorithm which will be used to analyze 
AIRS/AMSU/HSB data on the EOS Aqua platform. 
The following sections will describe the basic- 
methodology used to estimate cloud cleared AIRS 
radiances, which are subsequently used to retrieve 
surface and atmospheric geophysical parameters other 
than cloud parameters as well as to derive the effects ol 
clouds on the channel noise covariance matrix; describe 
the inversion methodology, which makes strong use of 
the channel noise covariance matrix and is applicable to 
solving for all the geophysical parameters including 
cloud parameters; and show sample results from AIRS 
Science Team simulations. 



2. CLOUD CLEARING METHODOLOGY 

Clouds have a significant effect on observed infra-red 
radiances, and can have smaller but non negligible 
effects on microwave observations as well. Therefore, 
an accurate treatment of the effects of clouds on the 
observed AIRS radiances is critical to obtaining 
accurate soundings. There are three basic approaches 
for treating cloud effects on the IR observations: look 
for clear spots and therefore avoid the problem, attempt 
to solve for the radiative effects of clouds directly in the 
inversion process; and attempt to infer what the 
radiances in the clear portions of the scene would be, 
called clear column radiances, from observations in a 
number of adjacent fields of view. An example of the 
first approach is given by Cuomo eLal. (1993). Eyre 
(1989a, 1990) has used the second approach in 
simulation by assuming an unknown homogeneous 
amount of black clouds at an unknown pressure, and 
attempted it with real TOVS data as well (Eyre, 1989b). 
Our approach, like that used in Susskind etaL (1997), is 
of the third type and is an extension of that used by 
Smith (1968), and Chahine (1974, 1977). The 
advantage of this approach is that it does not have the 
clear sky sampling bias of the first approach, nor does it 
require the ability to accurately model the spectral 
emissive, reflective, and transmissive properties of the 
clouds, and their dependence on the vertical 
microphysics and geometry, as required by the second 
approach. The key assumption made in the third 
approach is that while there may be many types of 
clouds in the different fields of view, the radiative 
properties of a given type of cloud are identical in all 
fields of view, which differ only in the relative amounts 
of these cloud types. Fields of view containing clouds 
with the same optical properties but at different heights, 
or clouds at the same height but with different optical 
properties, can be considered as having multiple cloud 
types. The other key assumption of this approach is 
that the fields of view have the same characteristics in 
the clear portions of their scenes, with unknown 
temperatures, humidities, etc. that we are trying to 
solve for. We have used analogous assumptions in 
analyzing 22 years of TOVS data on board the NOAA 
operational satellites (Susskind eLal, 1997) and shown 
that retrieval accuracy does not degrade appreciably 
with increasing cloud cover (Chahine and Susskind, 
1989). Analogous assumptions are made by 

NOAA/NESDIS in production of their clear column 
radiances used in generation of operational 

HIRS2/MSU retrievals (McMillin and Dean, 1982). 

Using these assumptions, Chahine (1977) has shown 
that in the case of K-l cloud formations, observations in 


K fields of view are needed to obtain channel i clear 
column radiances R i according to 

Ri =Ri 1 + I Ilk (Ri.1 -Ri,K+l-k) (1) 

k=l 

where R;.k is the channel i observation in field of view 
k. We have found it is advantageous (as suggested by 
L. McMillin) to extrapolate the radiances in the K fields 
of view according to a similar equation of the form 

R i = R iAVG + £ T lk( R i.AVG'' R i,k) 
k=l 

where Ri,AVG is f* ie average radiance of all fields of 
view. Optimal values of q k will give true values of R, 
up to instrumental noise effects. 

Cloud formations should be distinguished from cloud 
types. For example, if three fields of view are 
considered, and two cloud types exist, with cloud top 
pressures at 300 mb and 700 mb, and the respective 
cloud fractions as seen from above are (10%, 20%), 
(20%, 40%), and (30%, 60%) in each field of view, 
then only a single cloud formation exists with cloud 
fractions of 30%, 60%>, and 90% in each field of view 
respectively. If instead, the third field of view had 
cloud fractions of 30% and 65%, then 5% of a second 
cloud formation exists in the third field of view only. 
The above discussion applies only to cases in which the 
upper cloud type is opaque, and a portion of die scene, 
as observed from above, corresponds to cloud type 1, 
cloud type 2, or the surface. If the upper cloud type is 
semi-transparent, then a portion of the scene can 
correspond to cloud type 1 overlaying the surface, cloud 
type 1 overlaying cloud type 2, cloud type 2, and the 
surface. In such a case, three cloud formations will 
exist in general even if the relative amounts of each 
cloud type are as initially stated above. 

The methodology we use to determine % is general for 
handling up to K-l cloud formations. The simulations 
done by the AIRS Science Team, and shown in this 
paper, used essentially two cloud formations of gray 
clouds with differing amounts of clouds at 2 discrete 
levels in each of the 9 AIRS footprints within an AMSU 
A footprint. The cloud spectral emissivities and cloud 
top pressures were allowed to vary slightly between 
fields of view, however. Surface skin properties also 
had some variability between fields of view. This 
allows for multiple degrees of freedom within the 3x3 
array of AIRS spots in a single AMSU A footprint. 
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Susskind et at. (1998) used the 9 AIRS spots within an 
AMSU A footprint to construct 3 fields of view used to 
determine 2 values of q to be used in Equation 1. 
Field of view 1 was comprised of the average of the 
observations in the 3 wannest spots in an 8 pm 
window channel, and field of view 3 was the average of 
3 coldest spots. We now use all radiances in all spots 
separately and determine 9 values of q . Given q k * 
clear column radiances for all channels can be obtained 
from Equation 2. As in Susskind et al - (1998), we 
determine the values q and from observations in a 
selected set of I(= 76) cloud filtering channels which 
are primarily in between lines in the 15 pm ( 0, band 
and in the 4.2 pm C0 2 bandhead region, with some 
additional channels in the window regions. If, for each 
channel i, one substitutes an estimate of Ri.clr for Ri 
in Equation 2, this gives I equations for K unknowns. 
The unconstrained weighted least square solution to this 
multilinear problem is given by 


too high and incorrect values of q k would be obtained 
which would reconstruct too high values of R , . To 
avoid this, we make sure that the profile used to 
estimate Rj.cLR » s consistent with observations in all 

AMSU A and HSB channels, thus insuring an unbiased 
temperature and moisture profile over coarse layers in 
the atmosphere. It would be a mistake to use an 
analysis or a forecast field directly to compute R^clr 
because this field, while potentially accurate, could be 
biased in the vertical. 

The iterative methodology to determine clear column 
radiances consists of four passes to determine q“ (n = 
1, 2, 3, 4), using four sets of conditions, described later, 
to estimate R IVlr • * n which R i.ci.R and hence q , 
become increasingly more accurate for each iteration. 
Each set of conditions has its own N n , reflecting 
expected errors in R S|clr - Ri,i • The diagonal term of 
the noise covariance matrix is modeled according to 


q K;d = [AR 'N -1 Ar] kxK ARTT 1 AR CL r (3) 

where AR is a IxK matrix with ARj k =Ravg ~Ri,k' 
ARclr >s an Ixl matrix given by 
ARi,CLR =Ri,CLR- R i,AVG, and N is an Ixl channel 
noise covariance matrix. 


N-l = NEANj + 
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The 9 radiances Ri.k are observed at 3 different zenith 
angles. Having observations at different zenith angles 
will cause additional contributions to AR which are not 
due to differences in cloud cover. To remove these, we 
adjust all observed channel radiances to what they 
would have been if taken at die central zenith angle of 
the 3x3 array of AIRS spots according to Goldberg el 
aL 2002. From now on, AR refers to adjusted 
observed radiances. 


8R, 8q(P l 
aq(P) a 
. q 

and the off diagonal term is given by 


N,J 3T 3T V 3e„ 
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The key to the accurate determination of q is obtaining 
the best estimates of AR^clr* along with an accurate 

treatment of the noise covariance matrix N. As in 
Susskind etaL (1998), we assume the noise in channel i 
used to determine q is dominated by errors in 
ARi.clr- The values of ARj, C LR which we use to 
determine q (and R , ) are iterative and are computed 
based on the current best estimate of all relevant surface 
and atmospheric properties. For optimal results, it is 
important for the estimated geophysical parameters to 
be unbiased over large regions of the atmosphere. For 
example, if the estimated temperature profile were 
uniformly too warm, values of ARjclr would all be 


where NEANj is the channel i instrumental noise and 
the remaining terms are contributions to errors in the 
computed value Ri.clr resulting from errors in 
estimated surface skin temperature, surface spectral 
emissivity, surface spectral bi-directional reflectance of 
solar radiation, and temperature and moisture profile 
respectively. The partial derivatives are determined 
empirically by computing the radiance using the current 
estimate of each parameter and recomputing it after a 
small change in that parameter. In Susskind et_aL 
(1998), the uncertainties, such as 5T s n , are specified so 
as to be indicative of the expected errors for that 
parameter in pass n. We now predict these errors on a 
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profile by profile basis for each pass by propagation of 
expected sources of error through the retrieval process 
in a manner to be described later. A principal source of 
retrieval error arises from errors in the reconstructed 
clear column radiances. These errors propagate into 
degraded estimates of all the variables shown in 
Equation 4. 


One eigenvalue X k is always zero because only 8 
linearly independent values of AR^ exist. In 
transformed space. 


~ ^ max t 

R| = Ri.AVG + XCk AR;,k 

k=l 

and the solution for £ k is given by 


( 10 ) 


Selection of Optimal Fields of View 

The effects of instrumental noise on the clear column 
radiances will in general be amplified from single spot 
noise values because the clear column radiances are 
expressed as a linear combination of the observations in 
different fields of view. If there were no other sources 
of error, the diagonal term of the clear column radiance 
noise covariance matrix in a given pass in Equation 2 
would be 

[SR-SR'Jii-NEANf-ACTik) 2 (5) 

where A(q k ) is the noise amplification factor, given by 


9 ( 1 

9 > 



1+ Xtv 


k=lV 9 
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A (r| k ) is approximately equal to [X% 2 ]' /2 because 

the first term, containing the factor 1/9, is small. It is 
desirable to find an accurate expression for clear 
column radiance which minimizes a(tu). We can 
minimize A (f| k ) by expressing Equation 2 in terms of 

radiances in an optimal set of fields of view, given by 
linear combinations of die original set. The optimal 
A(t]k) can be found by transforming the original 

contrast fields, to a new set, AR k , according to 
AR^ k =XU k , k '-AR iik ' 0) 

k' 

where U is the unitary transformation which 
diagonalizes AR' N 1 AR 

ju'-(AR'-N -1 AR)u]k,k'=kk -8 k ,k'. (8) 

This is equivalent to having selected 

R k = RaVG “XUk.k'O^AVG ~Rk')* ^ 

k' 


£ k = A. k (AR T ' NT 1 ARclr) (ID 

t T 

where ARkk is the transpose of ARu 

It is apparent that large eigenvalues >- k imply low 
values of Ck while small eigenvalues imply large (and 
undesirable) values of Ck . The eigenvalues themselves 
indicate the degrees of freedom in the radiances in the 
different fields of view corresponding to the different 
number of cloud formations. Typical cloud formation 
eigenvalues are the order of 1000. We discard all 
eigenvalues less than 25 and set Kmax accordingly, 
with the constraint that K max is never greater than 4. 
We also do not include any eigenfunction whose 
eigenvalue is less than the uncertainty in C k . given 
later in Equation 13. This reduces the noise 
amplification factor by suppressing noise in the solution 
K 

of Ilk = 3T Ukk'-Ck'. resulting in lower values of 
k'=l 

bl- 
under certain pathological conditions, one or more 
cloud formations may not result in significant 
eigenvalues of ART4 ^ AR and cannot be solved for, 
resulting in a poor solution. The most obvious example 
of this is a single cloud formation with a constant cloud 
fraction in each field of view. Here AR is comprised 
of noise only. The most common examples of this are 
all fields of view are clear, which is a benign case, or 
all fields of view are overcast, which is a case which 
must be otherwise identified and rejected. Likewise, 
with two cloud formations, if the lower cloud deck is 
overcast, a proper resonstruction of the clear column 
radiances cannot be obtained. In this case, if the cloud 
fraction of the upper cloud in fields of view k is (*ik> 

then the lower cloud fraction as seen from above, ct2k » 
isl-auc- In general, if a 2 k = A + Baj k forallk.then 
cloud formation 2 will have a zero eigenvalue of 
ARTT'aR up to noise effects. The benign case occurs 
when A=0, corresponding to a truly single cloud 
formation. 
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Contribution of clouds to the retrieval channel noise 
covariance matrix M 

The basic retrieval methodology described in the next 
section requires a channel noise covariance matrix M 
representing channel correlated errors in the terms 
(Ri -Rf 1 ) and - Rj 1 ) where Rf 1 is the radiance 

computed for channel i based on the m tb iterative 
solution. The channel noise covariance matrix is the 
sum of two parts, resulting from noise in the 
reconstructed clear column radiances 5Ri with noise 
covariance M, and noise in the computed radiances 
5Rj m due to uncertainty in the parameters assumed 
known, with noise covariance M . 

M*j = [SR8R']ij is the expected noise covariance 

matrix for the channel clear column radiances. The 
noise in Rj obtained from Equation 2 has two parts, 
arising from instrumental noise NEANj, and from 
cloud clearing errors coming from errors in £k . Errors 
in will cause channel correlated clear column 

radiance errors. Clear column radiances for those 
channels affected by clouds will have this additional 
error due to errors in £ . For the AIRS instrument, the 
channel noise is spectrally uncorrelated, giving the Final 
result 

[SR8R']ii = NEANf A(nk) 2 + [AR T SC&;'AR r ] u 


[SRSR'Jjj =[AR T 8C8C'AR T ] lj (12b) 

where 8£8£' is the error covariance of £. If N, as 
defined in Equation 4, is indeed representative of the 
noise in the determination of rj , then it can be shown 
(see Equation 38) that 

[S£ 8£'] kk' = [ar T N - 1 AR T ]~*= A. k -1 8kk' ■ (13) 

In the special case for which we determine that channel 
i does not "see” the clouds (he., stratospheric sounding 
channels or tropospheric sounding channels peaking 
significantly above the highest cloud top), the clear 
column radiance is best described as the average 
radiance in all fields of view. For these channels, the 


scene appears to be clear and we can define effective 
values of q CLR for “clear” channels as q k CILR = 0 for 
all k. For these channels (see Equation 6), 

a(,F- r H (l4) 

which is a noise reducer. For “clear” channel i, one can 
write 

My =~NEANf 8,j (15) 

where j is any other channel and 8jj is the Kronecker 
delta function. 

For a channel to be determined not to see clouds, it 
must be included in a list showing a 95% probability of 
not seeing a cloud, which is pre-computed as a function 
of cloud top pressure and zenith angle. In addition, the 
standard deviation of the radiances in the 3x3 array of 
AIRS spots must be less than twice the channel noise. 

For channels which see clouds, the clear column noise 
covariance can now be expressed as 

M y = NEANj NEANj A(ri k ) 2 Sy + T (AR& ARj k X£ ! ) . 

(16) 

Errors in clear column radiances can be larger than 
predicted by Equation 16, however, because Xk is just 
an estimate of (8^8^)^, . Moreover, Equation 16 does 

not take into account contributions to the noise 
covariance matrix arising from higher components of £ 
not solved for (k > K^) as well as fitting errors due to a 
poor first guess. Another estimate of the error in the £ 
parameters can be computed using weighted radiance 
residuals in the channels used in the cloud clearing 
retrieval, Rj,cLR“^i If we take Ri,CLR - Ri as toe 
uncertainty of AR^clr* toen using Equation 11, we 
estimate the uncertainty in £ k according to 

^(aR^Nu'^R^clr-R 1^(17) 

which we evaluate for all significant functions k with 

~3 

^k > 10 . This includes eigen functions with 
< 25 and therefore not included in the solution for Rj . 
For values of k < K max , we take 

[5C8CV=MAx[v^SC'] kk ] (18) 

and for values of k between K niui and K ik (significant 
eigenvalues > 1 0 1 ) we set 
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[§C5C']kk=[6C6C']kk 


( 19 ) 


and write 

My = NFAN, NEANj A(q k ) 2 
K . . (20) 

+ I AR,‘ k AR^ \KK] kk 
k = l 

One can think of Equation 20 in terms of a different 
effective noise amplification factor Aj e ff l° r eac ^ 

channel i 

M a = NEANj 2 A 2 i,eff (2D 


where 


Aj,eff - 




( 22 ) 


The channel effective noise amplification factor is 
largest for channels which see the surface and have 
potentially large values of the scene contrast AR lk 
We find it convenient to define an effective noise 
amplification factor relevant to the surface channel 
retrieval step as the RMS value of A , e ff over all 
NSURF infrared channels used in the surface retrieval 


step 


A eff = 


1 

NSURF 


NSURF 

I 

i=l 


it 12 


Ai.eff 


( 23 ) 


step and using previously determined variables as fixed 
but with an appropriate uncertainty attached to them 
which is accounted for in the channel noise covariance 
matrix M. The objective in each step is to find 
solutions which best match the observations for a select 
set of channels, bearing in mind the channel noise 
estimates. The "observations" in steps a-d are the clear 
column radiances as determined from Equations 10 and 
1 1, with values of C CLR used for appropriate channels. 
The cloud parameters determined from step e are found 
so as to be most consistent with the actual observed 
radiances and the clear sky geophysical parameters 
determined from steps a-d. Steps a-d are ordered so as 
to allow for selection of channels in each step which are 
primarily sensitive to variables to be determined in that 
step or determined in a previous step, and relatively 
insensitive to other parameters. Separation of the 
problem in this manner also allows for the problem in 
each step to be made as linear as possible. Steps a-e 
are all solved for in a completely analogous manner, 
linearizing the problem about initial guess parameters 
and iterating the solution until convergence is reached. 
In general, these linear equations are ill conditioned and 
require some form of stabilization, which is commonly 
based on an estimate of the accuracy of the a-priori 
information obtained in the first guess or background 
field (Rodgers, 1976; Hanel et _al, 1992). The 

methodology we have developed, described in the next 
section, relies exclusively on the signal to noise of the 
observations to indicate the degree to which the 
information contained in the radiances should be 
believed, and does not involve use of an estimate of the 
accuracy of the background field. 


Very large values of A eff can arise when 5£8£' is large 
(A ef f is sometimes 100 or more) and indicate a large 
uncertainty in the determination of the clear column 
radiances. These large uncertainties are sometimes 
caused by hidden, or nearly hidden cloud formations, 
and often correlate with poor solutions. 

3. BASIC RETRIEVAL METHODOLOGY 

The basic retrieval methodology is the same as that of 
Susskind et al. (1998) and is reviewed below. After a 
start up procedure to determine the clear column 
radiances, we use A1RS/AMSU/HSB data to retrieve: 

a) surface skin temperature, surface spectral emissivity 
and surface bi-directional reflectance of solar radiation; 

b) atmospheric temperature profile; c) atmospheric 
moisture profile; d) atmospheric ozone profile; and 

e) cloud properties. These steps are done sequentially, 
solving only for the variables to be determined in each 


Iterative least squares solution to the non-linear 
problem 

The solution to each of the five steps described above is 
done in the form 

X m = X m-t + !F f AA? = X° + IF, A? (24) 
r=t i = i 

where X m is the m th iterative state, F f is a set of L 
functions, and 

Af = Af~ l +AAf. (25) 

A A™ is determined each iteration so as to minimize the 

residuals A0 m , weighted inversely with respect to 
expected noise levels, for die channels used to 
determine Af . The residual for channel I is defined as 
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where Rj is the reconstructed clear column radiance, 
r 111 is the radiance computed from the m iterative 
parameters, and 0™ is the brightness temperature 
corresponding to Rj”. The m' h iterative residual for 

channel i is attributed to errors in the coefficients 5A m 
and to the noise effects, 

Ae‘ n =lS i ,5Af + 0 1 , ( 27 > 

e 

where Stf is an element of the sensitivity matrix, or 
Jacobian, given by 


profile, we are holding fixed surface parameters, 
moisture profile, and ozone profile. We currently write 



and 


3R. 
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(28) 
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(30b) 


and the noise factor 0, for a given case has two parts: 
errors in observed clear column radiance 8©j which are 
affected by instrumental noise and cloud clearing 
errors, and computational noise 80; . 

In the simulations done thus far, we have assumed 
perfect knowledge of physics, i.e., if we know all of the 
variables exactly, we can compute exact noise free 
radiances. Nevertheless, the transmittances depend on 
the variables to be solved for. Therefore, 
computational noise exists. Computational noise, 
arising from errors such as a low (high) estimate of 
atmospheric water vapor, will produce noise that is 
spectrally correlated. Instrumental noise is spectrally 
uncorrelated but cloud cleared radiance errors are 
correlated. Each retrieval step in pass n uses an 
appropriate noise covariance matrix 



where M was defined in Equation 20. 

The matrix M in Equation 29 represents channel 
correlated uncertainties in the computed radiances R , 
and Rj" based primarily on uncertainties in the 

parameters being held fixed in a given retrieval step. 
For example, when we are solving for temperature 


The terms in Equations 30a and 30b do not depend to 
first order on the variables being solved for or the 
iterative state m. The term 0.1 2 included in Equation 
30a is taken to represent additional uncertainties in 
computed brightness temperatures based on the 
imperfect knowledge of the variables being solved for, 
as well as potential spectroscopic errors. The 
methodology used to predict and propagate errors such 
as 8T S for use in the computation of M will be 
discussed later. These terms are analogous to the terms 
in Equation 4, but uncertainty in 0 3 profile was not 
included in Equation 4 as it did not prove to have a 
significant effect on the solution and the calculation is 
computationally expensive. 

A pplication of a constraint H 

The standard constrained solution (Hanel etjU-. 1992) 
to this problem is given by 

AA m =[s m 'M"'S m +H m ] 1 S m 'M- 1 A0 m (J1) 

= D m A0 m 

where H m is a stabilization matrix. Without 
stabilization, SAA would minimize the weighted 

residuals A0'M _1 A0, but die matrix elements of D 
might be large. This is undesirable as it amplifies errors 
in A0 m in determining A A m . The key to optimization 
of the solution lies in accurate treatment of the tenns 
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A0 , S and M; a judicious choice of the functions F and 
channels i; and optimal treatment of the constraint 
matrix H. Hanel et al. (1992) and Rodgers (1976) have 
reviewed several methods of constraining the ill - 
conditioned inverse problem. In the minimum variance 
approach (Rodgers, 1976), H is taken to be the inverse 
of the a priori error covariance. If the statistics of both 
the measurement and a priori are Gaussian, the 
maximum likelihood solution is obtained. If die a 
priori covariance is taken to be H = yl . the maximum 
entropy solution is obtained. Other forms of H include 
the first or second derivative formulations (Twomey, 
1963) that force a smoothness constraint on the 
solution. The solution can also be constrained by the 
relaxation method (Chahine, 1968) and by the Backus 
and Gilbert (1970) method. 

The minimum variance and maximum likelihood 
solutions are often considered to be “optimal.” 
However, if the a priori error covariance is not known 
or estimated incorrectly, the solution will be sub- 
optimal. If the a priori errors are underestimated, the 
solution could be overconstrained. This could 
potentially create biases in the retrievals. The biases 
may mask small trends in die retrieved data that one 
may be trying to extract. The approach described here 
attempts to keep the effects of instrument noise at a 
tolerable level without assumptions regarding the a 
priori data error covariance. 

Our objective is to determine a constraint matrix that 
affects only the pieces of information not well 
determined by the radiances. This involves use of an 
optimal set of functions G, related to F by a unitary 
transformation 

G = FU ( 32 > 


eigenmode G, . The unconstrained solution, widi both 
Hy = 0 and with no background correction (80 n = 0j , 
is given by 

AB™ (0)= (*? J 1 (u m S m 'M- 1 ) f A0 m 

= (^ l ? J 1 d™A0 m 05) 

where d" 1 is the f th row of U ra 's m 'M -1 . Low values 
Of indicating ill-conditioned variables would 
lead to large coefficients of A0 n in the determination 
of AB™ if the solution for the coefficients of these 
functions were unconstrained. We therefore only 
constrain the solution of those functions with low 
eigenvalues and set H™' = AA^ 5^' . The constrained 

solution is now given by 

t m 

4Brte>T^ X -^4B!’(0)=0?4B?(0) (36) 

^ ' A ^ + AA^ 

where O is a damping factor equal to unity if no 
constraint is applied and zero if AXg is infinite. If 
d>/ =1 , all the information about G , obtained from the 
radiances is believed. Otherwise only Of of the 
information in the radiances is believed, and (l-Of) 
of the a-priori information is believed. 

The objective of damping a mode is to reduce 
propagation of noise 8B/(AA./ ) which is given by 

8B? + A*?)"' (u m S m 'M“ l )*S> m (37) 

where 80 m is the noise in A©i . A statistical estimate 
of 8B™ over an ensemble of profiles is given by 


in terms of which we could write 

X m = X m_l + C-AB m = X" 1-1 + FIJAB” = X m_1 + FAA m . 

(33) 

In terms of the functions G, Equation 31 becomes 
AB m = U'AA m = (U m/ S m 'M“ l S m U m + H m ) 

• u m S m M -1 (A0 m -80 m_1 ) . (34) 

A new term 80 m_1 has been included in Equation 34 
which is a background correction term that is zero in 
the first iteration and will be discussed in detail later. 
The optimal transformation matrix U m is chosen so 
that u m S m M -1 S m U m is diagonal with eigenvalues 
X ( . The inverse of each eigenvalue is the variance of 


SB, (A x ) = [® m SB' m ] « 2 
[u^S™ M- , 5080'M- 1 S m U m ] ^ 

= _£ii — = 

V? + A A? V J 


1/2 
I ee 


(38) 


to the extent that M accurately represents 8080' , the 
noise covariance of A0 . We assign a noise 
propagation threshold AB max for each type of retrieval 
(surface properties, temperature profile, etc.) and set 
AXf=0 if 8Bf_(0)<AB max . Otherwise, we obtain 

A\( such that 8Bf(AX,) = AB max * 
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Formulation of the background term 


The need for an iterative process arises because the 
radiative transfer equation is not linear. In every 
iteration, we recompute @i , as well as 
S m , U m and X m . If the solutions were completely 
linear, and we applied no damping, then 

A0 m+1 ( O ) = 0 - Q m+1 (0)= A0 m - S m U m AB m (0) 

(39) 

y 

and AB m+l (0) would be determined to be zero 
because AB m (0) would have already minimized the 

residuals A0 m . 

The residual A0 m+1 is not zero however, both 
because © + (0) is not given exactly by 
e m + S m U ra AB m (0) as a result of non-linearity, and 
because AB? * AB? (0). As a result of applying 
AB™ rather than AB? (0), we obtain 

A0 m+1 « A0 m+1 (O) + S m U m [AB m (0) - AB m ] 

= A0 m+1 (0)+ 50 m < 40 > 

In Equation 40, A0™ +1 (O) represents the portion of 
A0 m+1 that is due to effects of non-linearity on the 
solution, while 80 m represents the residual portion 
of A0 m+1 due to the effects of damping in iteration 
m. The second term is zero for undamped modes and 
increases in significance with increased damping. 
This term is also zero for all modes in the first 
iteration. We only want to include the effects of non- 
linearity on A0 m in the iterative procedure used in 
the determination of A® . Therefore, the 
background term to be used in Equation 34 is given 
by 

g0m -l =s m-l ym-1 ( 0 )_ AB m_1 ] (41) 

and we solve for AB? according to 


ab?=(X7+a x?)' 1 U m ' s m ' m- 1 

. [A0 m -80 m_1 ]=O m AB?(Q )- (XT + AA-T)' 

. [u m 'S m ' M" 1 S 1 "' 1 U m_1 (ABr'(0)- AB^T 1 )] 

(42) 

where AB? 1-1 is the value of AB f which was applied 
in iteration m-1. Inclusion of the background term in 
Equation 42 insures second order convergence along 
the lines discussed by Rodgers (1976) with regard to 
treatment of the a-priori term. 


Convergence Criteria 


In solving Equation 42, we are attempting to find 
solutions to the radiative transfer equations which 
minimize weighted residuals of observed and 
computed brightness temperatures, corrected for the 
background term. To test convergence of the 
solution, one should monitor the weighted residual 


R = [(A0- 80)’ V' V (A0 - 50) 


, (43) 


where the weight matrix V accounts for noise effects 
on the channel residuals, as well as the relative 
information content of the channels with regard to the 
variables being solved for. For example, if a channel 
(or linear combination of channels) carries little 
information content in terms of signal to noise, it 
should be given little weight in the estimation of the 
residual in Equation 43. An appropriate choice of V, 
expressing the information content of the channels, 

would therefore be V = (X* + Ak f ) (u S24 ), in 


which case we obtain 


r=[ab'ab] 1/2 . 


(44) 


As shown in Equation 44, a reasonable way to 
determine if the solution has converged, in terms of 
weighted residuals of observed minus computed 
brightness temperatures, is to see if the solution has 
converged in terms of the iterative changes in the 
solution itself. Initially, we set ABj =0 if 0^<0.05, 

that is, coefficients of very heavily damped components 
with little information content are not believed at all in 
any iteration. The solution is said to have converged 
when the RSS value of AB? is less than 10% of the 
RSS value of 5B m for all components not set equal to 
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zero. The iterative procedure is also terminated if the 
RSS value of AB™ is not less than 75% of that of 
AB™ -1 for die non-zero components. This indicates 
the solution is not converging rapidly enough and may 
be responding primarily to unmodeled noise. The 
iterative procedure, which typically converges by 3 
iterations, is carried out analogously for all retrieval 
steps. 


Variables and functions for retrieval steps 

As shown in Equation 24, all steps involve expression 
of the basic variables X in terms of a set of functions F. 
In the temperature profile retrieval step, the temperature 
perturbation functions F(P) are 24 trapezoids piecewise 
linear in log of pressure, spanning the pressure range 
.016 mb to the surface, with a perturbation of 0.5K 
between pressures P< and P/ _j , and 0° at 
P f +i and Pp- 2 . In the top and bottom functions, the 
top or bottom portion of the trapezoid is missing. The 
Jacobian s!? is obtained numerically by computing the 
channel i brightness temperature for the m th iterative 
temperature profile T m (P) and subtracting it from the 
brightness temperature computed with everything else 
fixed but perturbing T m (P) by one unit of F ( (P) . 
With regard to water vapor and ozone profiles, we 
express solutions in the form 

q m+, (P)=q m (P)[l + lA,F,(P)] (45) 

with the functions (L = 11 for water vapor and 8 for 
ozone) and methodology for computation of the 
Jacobian being completely analogous to those for 
temperature profile. In the case of surface variables, 
the functions are a discrete value of surface skin 
temperature, as well as 9 triangle functions in the 
frequency domain dealing with perturbations of surface 
emissivity and 3 with surface bi-directional reflectance 
of solar radiation. The total precipitable water also can 
be adjusted by using Equation 45 with a single function 
which is constant as a function of height. The window 
channels are sensitive to boundary layer water vapor 
but not higher level water vapor. Adjustment of total 
precipitable water is used in an intermediate retrieval 
step, done before the water vapor retrieval step using 
AIRS channels, and improves total low level water 
vapor at the expense of upper level water vapor. Table 
1 shows the pressure levels used in the temperature 
profile, moisture profile, and ozone profile retrieval 
steps. 


Table 1. Trapezoid Function Endpoints (mb) 


Temperature 

Moisture 

Ozone 

Retrieval 

Retrieval 

Retrieval 

0.016 

0.016 

0.016 

0.714 

170.1 

20.92 

1.297 

272.9 

51.53 

2.701 

314.1 

71.54 

4.077 

343.6 

103.0 

8.165 

407.5 

142.4 

16.43 

496.6 

300.0 

23.45 

617.5 

surface 

39.26 

706.6 


56.13 

852.8 


71.54 

surface 



96.11 

125.6 

160.5 
212.0 
272.9 

343.6 

424.5 

496.7 
596.3 
661.2 

753.6 

878.6 
surface 


Selection of channels 

While AIRS has 2378 channels, it is neither necessary 
nor optimal to use all the channels in the retrieval 
process as the information content of these channels is 
highly redundant. Therefore, computational time can 
be lowered by limiting the number of channels used. In 
a given step, it is preferable to use channels which are 
primarily sensitive to the variables being solved for, 
while relatively insensitive to variables not yet 
solved for. We also find it desirable to use channels 
with sharp localized weighting functions. Kaplan et al. 
(1977) show that channels with sharpest weighting 
functions lie either in between absorption lines or on 
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the band head of the 4.3 |iin C0 2 band between 2378 
cm 1 and 2390 cm 1 . The first set of channels have sharp 
weighting functions because of a rapid increase of 
absorption coefficient with increasing pressure, while 
the second benefit from a rapid increase of absorption 
coefficient with increasing temperature in the 
troposphere. Such channels form the basic set used for 
temperature sounding. Channels in between water 
vapor absorption lines also produce very sharp 
weighting functions which are preferable for water 
vapor sounding and also useful for temperature 
sounding if the water vapor distribution is known 
accurately. Channels between absorption features are 
by definition less opaque than nearby channels situated 
on absorption features, and may not have sufficient 
opacity to be sensitive to either temperature or 
constituents at high enough levels in the atmosphere. 
For temperature profile, we select channels in the 
C0 2 Q branch at 667 cm 1 , which do not have sharp 
weighting functions but are sensitive to temperature 
variations up to 1 mb. We do not select channels in the 
most opaque portion of the 4.3 pm C0 2 band because 
these channels are sensitive to effects of non-local 
thermodynamic equilibrium. For water vapor, we 
selected a few channels on the peaks of some of the 
strongest absorption features in the 6.7 pm water vapor 
band to increase the sensitivity to stratospheric and 
upper tropospheric water vapor. Window channels are 
highly redundant with each other and have been 
selected generally on and off closely lying weak 
absorption features in the spectral regions from 755 
cm ' - 980 cm 1 , 1070 cm 1 - 1240 cm 1 , 2180 cm 1 - 
2192 cm 1 , and 2390 cm 1 to 2665 cm 1 . Cloud filtering 
channels are generally a subset of the temperature 
sounding channels which are sensitive to the 
troposphere. Our sounding methodology involves two 
temperature profile retrieval steps, one (temp 1) before 
the water vapor retrieval step, and the other (temp 2) 
subsequent to it. In temp 2, we include a number of 
channels in the water vapor absorption band which 
produce sharp temperature weighting functions. These 
channels are treated as “noisy” in the channel noise 
covariance matrix to the extent that the prediced 
uncertainty in water vapor distribution produces an 
appropriate uncertainty in their computed brightness 
temperatures. The location of all channels used are 
shown in Figure 1. We use 53 channels in the surface 
temperature retrieval, 147 channels in the first 
temperature profile retrieval, an additional 7 channels 
in the second temperature profile retrieval, 66 channels 
in the water vapor profile retrieval, and 23 channels in 
the ozone profile retrieval. Some channels are used for 
more than one purpose. Channels also exist which can 


be used for retrievals of profiles of CH 4 , CO, and CO,. 
These will be described in a future publication. 

Table 2 shows the eigenvalues and damping factors for 
the second pass temperature profile retrieval, the water 
vapor retrieval, and the ozone profile retrieval for a 
typical case. Coefficients of eight temperature profile 
functions are undamped, and those of two more 
functions are only slightly damped, giving about 9 
pieces of information about the temperature profile 
being contained in the radiances. Roughly 4 1/2 pieces 
of information about water vapor are contained in the 
radiances, and roughly 1 1/2 pieces of information are 
contained about the ozone profile in this case. 

Table 2. Sample Eigenvalues and Damping Factors 
(O > 0.05) 

Temperature Profile Water Vapor Profile Ozone Profile 



= 0.75 

AB max 

= 1.0 

AB max 

= 0.75 

X 

O 

X 

O 

X 

d> 

94.40 

1.0 

69.78 

1.0 

41.59 

1.0 

65.06 

1.0 

10.31 

1.0 

0.745 

0.419 

36.81 

1.0 

2.44 

1.0 

0.403 

0.278 

26.55 

1.0 

1.26 

1.0 

0.139 

0.079 

12.19 

1.0 

0.68 

0.668 



6.96 

1.0 

0.27 

0.270 



3.84 

1.0 

0.09 

0.095 



1.84 

1.0 

0.07 

0.069 



1.21 

0.685 





0.80 

0.447 





0.50 

0.281 





0.35 

0.1% 





0.29 

0.161 





0.20 

0.114 





0.11 

0.006 






Cloud parameter retrievals 

In performing cloud parameter retrievals, all other 
variables are assumed known within their estimated 
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Figure 1. AIRS spectrum showing channels used in 
different retrieval steps. Temperature sounding 
channels are red, ozone are green, water vapor are 
blue, and surface channels are orange. 

errors, allowing us to compute R^clr * The channels 
used are the subset of cloud clearing channels that are 
not sensitive to solar radiation reflected off the 
clouds. The cloud parameter retrieval algorithm is 
analogous to that of the other steps but slightly 
different. At this time, the cloud retrieval algorithm 
has been tested only for the case of assumed cloud 
spectral properties in order to determine cloud 
fractions and cloud top pressures for up to two layers 
of clouds. The method is easily generalizable to 
include cloud spectral emissivity by inclusion of an 
appropriate set of spectral emissivity functions as 
done in the surface parameter retrievals. With known 
spectral properties, cloud radiances Rj(P c ) can be 

calculated based on the surface skin 
temperature and atmospheric temperature-moisture- 
ozone profile, which have been retrieved from the clear 
column radiances and are "known”, as a function of 
unknown cloud top pressure P c . For two cloud layers 
(the method works for any number of cloud layers) we 
can write 

R ik = ( 1 ~ a lk “ a 2k) R i,CLR 

(46) 

+ CtikRi(Pcl )+ tX2k R i(Pc2) 


where Rik is the radiance computed for channel i in 
field of view k covered by (as seen from above) oqk 
fractional coverage of a cloud at P cJ and a 2 k °f a 
cloud at P c2 . In the above equation, we have assumed 
two types of clouds in each of the Fields of view k=l,9, 
with different cloud fractions in each field of view. All 
clouds were assumed to have a constant spectral 


emissivity of 0.9. In order to determine the variables 
Pci,Pc 2 ,an,ai 2 ..., we use observations in the 9 

fields of view for the subset of channels used to 
determine r\ which are unaffected by solar radiation. 
The noise covariance matrix N used to retrieve cloud 
parameters is identical to that used in Equation 4 to 
determine t| , but for the appropriate subset of channels. 


„ th , , mm 

Given the m guess cloud parameters aik . a2k . Pci . 
and Pc™ , we define 

Y,k =Ri.k-R|k=(R.,k-Ri,cLR) 

+ I^( R .,clr- R .( P cT)) (47) 

and obtain the iterative equation 


Y r 1 -Y i |T= I [(Ri,CLR 
j=l,2 LV 



Aa 


jk 



3P C , 



where the terms in the square brackets are the 
appropriate Jacobians, which are computed 
empirically as are all other Jacobians. It should be 

noted that if a j k (for all k) and/or — 1 (for all i) are 

Or c j 

small for a given P c j , that cloud top pressure will be 


contained primarily in a heavily damped mode and 
not be changed significantly from the initial guess. 


Error Propagation and Channel Noise Covariance 
Matrix 


Equations 4 and 30 contain terms such as 5T(P), 
indicative of expected errors in state parameters used in 
a given pass and step. These errors are case dependent 
and can be estimated by propagating expected errors 
through the retrieval system. In any iteration, the 
estimate of a parameter, such as T(P) m , is given by 
L 

T(P) m =T(P)f + Z F..A" 1 

J 1 r = t 1 (49) 

= T(P)j + (FU B m ). j 

where j is a discrete pressure level. There are three 
contributions to the expected error SlXP)™. The first 
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contribution comes from the null space error, arising 
from the error of the first guess in the space outside that 
of the L functions used to expand the solution. The 
second component arises from errors in the coefficients 
B ra . The last contribution arises from the damping of 
the solution in which (1-0 ) of the first guess (or 
previous iteration) is believed for each eigen function 
G=FU. 


next retrieval step and used in Equation 51 to give 
5Ak° which is in turn used in Equation 52 to generate 
the uncertainty in parameter X for use in subsequent 

steps. 

For moisture and ozone profile, the form of the 
expansion is slightly different (see Equation 45) and 
we write 


Equations 4 and 30 contain the square of the expected 
error in state parameter X m , 8 X" 1 , which can be 

expressed in terms of errors in the expansion 
coefficients A according to 


8X"' 2 = 8X^ 2 + 1 F^SAjP 2 (50) 

k 

where 8X^ is the null space error and 8A m is the error 

in the coefficients A m used to represent X m . Errors in 
A arise both from errors in the B coefficients and 
errors in the damped portion of the m-1 iterative guess. 
In every step in the retrieval process, we begin with 

parameters X° having an uncertainty 8Xj°. The 

uncertainty of the microwave product first guess is 
specified based on expected errors, as is the null space 

error. Given 8X° , 8A° can be solved for according 
to 


8A k ° =[(F ,2 F 2 r 1 F ,2 (8X°-8X N ^j 
= |(F' 2 F 2 J l F' 2 (5X°)fj 


(51) 


In a given iteration, we can express 5A™ according to 




/ \ 

2 

r \ 

2 

SA? = 


« m 

rt,\ 

+ £ f 

U W (1 -♦ < -)ZU j< 8A“-' 

, i J 

R N 2 

+§A^ 


(52) 


where 



represents the predicted error in B™ 


due 


to noise propagation, and die second term represents the 
portions of the errors 5B^ m_1 of the previous iterative 
profile which are believed in the current iteration. 


Given 5 A™ from Equation 52 for die final iterative 

step, we compute the square of the corresponding 
profile error to be used in Equations 4 and 30 
according to Equation 50. This term is carried to the 


{ 5q m (P)] 

2 

f8q N (P)f 

- 

JC 5 




o J 


+IF k 2 (P)8A|T 2 

k 


(53) 


Surface spectral emissivity and bi-directional 
reflectance are analogous to temperature profile, as is 
skin temperature, in which case F is a number. The 
liquid water profile comes from the microwave product 
and is not iterated. We assume an error estimate of 
20% of the liquid water profile. In addition, if the total 
liquid water is less than 0.01 g/cm 2 , we consider the 
possibility that liquid water may have been missed due 
to an error in the water vapor microwave solution. For 
these low liquid water solutions, an alternative error 
estimate of (2*RH-l)*0.05*q, where RH is the relative 
humidity and q is the layer water vapor in mg/cm 2 , is 
considered and used if it is larger than 20% of the liquid 
water. The null space temperature error is taken as 
0.1K in the lower and upper atmosphere, increasing to 
0.2K near the tropopause. The null space error in 
percent is taken as 5% for water vapor and 10% for 
ozone respectively. 


Equation 52 is case dependent through the parameters 
<t> ( and which depend both on the S matrix, and 

more significantly on the M matrix. M contains 
contributions from clouds, M , and parameter 
uncertainty M . The uncertainties determined from 
Equations 52, 50, and 53 in turn are used in the 
computation of M (Equation 30) and N (Equation 4). 


Equations 50 and 53 give the magnitude of the 
estimated error in each parameter but contain no 
information about sign. If we assume all 8X(P) are of 
the same sign, we would overestimate the effect of the 
uncertainty on that parameter on the computed 
radiances. Bearing this in mind, when the derivatives 
in Equations 4 and 30 are computed numerically , we 
write 


-^-dX(P)= R(X(P) + AX(P))-R(X(P)) 
oX(P) 


(54) 
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where AX(P) is constructed by multiplying 3X(P) by a 
sine wave with a full period of six temperature profile 
functions in the case of uncertainty of temperature 
profile to be used in the humidity and ozone profile 
retrievals, and six humidity profile functions in die case 
of water vapor uncertainty to be used in the temperature 
and ozone profile retrievals. In the case of ozone 
profile, with only seven functions, we simply multiply 
the predicted uncertainty by 0.5. We have also found 
that in constructing die noise covariance terms in 
Equation 4, it is advantageous to set AX =03 8(X) for 
all profile terms. For surface parameters we take 
AX = 3X , as for the liquid water profile. 

Steps in the processing system 

The processing system used in this paper is comprised 
of a number of sequendal steps listed below. All steps 
start from the conditions found in the previous step, 
with appropriate computed uncertainty estimates, 8X , 
unless otherwise noted. 

1. Use as a starting point the microwave product 
which agrees with the AMSU A, HSB 
radiances (Rosenkranz, 2000). We follow this 
by a temperature profile retrieval using AMSU 
A radiances as well as AIRS radiances for 
channels diat never see clouds, followed by a 
water vapor retrieval using HSB channels and 
some AMSU A window channels. As part of 
the temperature profile retrieval, we also 
update the surface skin temperature and 
microwave spectral emissivity. 

2. Determine an initial from Equation 11 
using the atmospheric and surface parameters 
obtained in Step 1. We also perform a cloud 
parameter retrieval to help determine which IR 
channels are not affected by clouds. R{ is 

obtained using in Equation 10. 

3. Determine the first guess IR surface 
parameters and temperature-moisture-ozone 
profile using Ri based on a regression step 
using most AIRS channels (Goldberg et al., 
2002). Under some difficult cloud conditions, 
tins first guess is modified in a manner 
described later. 

4. Produce an improved temperature profile and 
microwave spectral emissivity starting from 


the surface and atmospheric parameters 
determined in step 3 using the AMSU A 
channel radiances and AIRS channel radiances 
which do not see clouds. The surface skin 
temperature is not updated as it is estimated 
better from AIRS radiances than can be 
determined from AMSU radiances. This is 
followed by an improved water vapor profile 
using HSB radiances. 

5. Determine taking advantage of the refined 
parameters. Also determine cloud parameters 
to decide which channels do not see clouds so 
as to average radiances in these channels when 
producing R 2 . R 2 is considerably more 
accurate than r{ because the surface and 
atmospheric parameters obtained from the 
AIRS regression step are more accurate than 
those from the microwave first product, 
especially the infra-red surface spectral 
properties which are not determined from the 
microwave retrieval. 

6. Perform a surface parameter retrieval using 
AIRS surface sounding channels shown in 
Figure 1, and AMSU A and HSB window 
channels. This produces a new skin 
temperature, IR and microwave spectral 
emissivity, and IR spectral bi-directional 
inflectance. It also includes adjustment of the 
entire water vapor profile by a single 
trapezoidal function which is constant in the 
troposphere and lower stratosphere. 

7. Determine £ 3 taking advantages of the refined 
surface parameters, and produce Rj and new 
estimates of cloud parameters. 

8. -11. Use R 3 to sequentially determine surface 

parameters, temperature profile, humidity 
profile, and ozone profile. These are called the 
first pass retrieved products. 

12. Update the temperature profile, using only 
AMSU A radiances and AIRS channel 
radiances insensitive to clouds. This profile is 
also used in the rejection criteria and is 
referred to as the test microwave only retrieval. 

13. Using the first pass retrieved products and 
updated temperature profile, determine Q , 
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final cloud parameters, and the final clear 
column radiances R? , which is a product of 
the system. 

14. Repeat steps 8 and 9 using Rj* to obtain the 
final product surface parameters and 
temperature profile. The initial guess used in 
the second pass surface parameter and 
temperature profile retrievals is identical to 
that of the first pass but all other parameters 
are updated, such as the clear column 
radiances, moisture profile, etc. The noise 
covariance matrix is also updated to account 
for better estimates of the other parameters. In 
addition, channels in the water vapor band 
which are highly sensitive to lower 
tropospheric water vapor are included in the 
final temperature profile step (but not the first 
pass) because an accurate moisture profile has 
now been retrieved. The moisture profile and 
ozone profile retrieval steps are not repeated, 
as no appreciable improvement in parameters 
resulted from further retrieval steps. 

15. Test solution for acceptance. If rejected, 
return to the AMSU/HSB retrieval starting 
from the initial guess, including AIRS 
channels insensitive to clouds, as the “final 
microwave only” product. Cloud parameters 
for rejected cases are based on this solution 
and were determined in step 2. 

16. Determine OLR and clear sky OLR using the 
appropriate solution for either accepted or 
rejected cases. 


Adjustment of the First Guess 

The first guess temperature profile T°(P) used in step 
4 is usually the result of the AIRS regression done in 
step 3, T reg (P) . Under most conditions, this is 
considerably more accurate than the microwave product 
T (P) used in step 1. However, under some difficult 
cloud conditions, a very poor regression can be 
obtained. In general, the regression temperature profile 
will degrade below 300 mb with increasing values of 
the effective noise amplification factor, and can be 
considerably poorer than the microwave retrieval, 
especially near the surface. The problem is 
compounded, in cases of large effective noise 
amplification factor, because temperature sounding 
channels sensing the lower troposphere will be treated 


as noisy due to a large contribution of the second term 
in the noise covariance matrix (see Equation 20). 
Consequently, eigen functions having high vertical 
resolution in the lower troposphere will be heavily 
damped and the poor vertical structure in the lower 
troposphere of the regression guess will be heavily 
believed. To help alleviate this problem, we construct a 
first guess temperature profile which is a linear 
combination of the T reg (P) and T M (P) below 300 mb. 

T°(P)=T reg (P) + A(P)(r M (P)-T reg (P)) (55) 

where A(P s ) = 0 if A e f f <2.5 and A(P S ) = 1 if 
A ef f >6.5, with values of A(P S ) linearly interpolated 
between 0 and 1 for intermediate values of A eff . In 
addition, A(P) = 1 for P<300mb, and A(P) is 
linearly interpolated in ^nP for intermediate pressure 
values between 300 mb and the surface pressure P s . 


Computation of OLR 

Outgoing Longwave Radiation (OLR) is computed 
from the AIRS products in a manner analogous to 
that used to compute OLR from TO VS (Mehta and 
Susskind, 1999a, 1999b). 


F = (l -ax -a 2 )F CLR + «t F cld (P Ci ) + a 2 F CLD (P C:! ) 

(56) 

where F<xr > the clear sky OLR, is the sum of 
contributions from 14 spectral bands each with 
effective surface emissivity Cj 
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f clr = K £ 
M 


e i B Vi (T s )x i (P s ) + | B Vi (t)— M ft»P 
ML QK[U 


= K^Fi.CLR 
i 


(57) 


and the band transmittances Xj(P) are computed at 
effective zenith angles 0* . The small term related to 
downwelling thermal radiation reflected off the 
surface and transmitted to space is neglected. 
Fcld(Pc) * s computed in an analogous way, in terms 
of the cloud spectral emissivity £i(P c ), anti 
assuming a cloud transmissivity of (l- £j(P c )) 
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£i (p c )B v . (T(p c mj(p c > 
dx 
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+( 1_e >^i>CLR 


(58) 
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The band transmittance s (P) are parameterized as a 

function of temperature, moisture, and ozone profile 
(Mehta and Susskind, 1999b), The spectral cloud 
emissivity was assumed to have a constant value of 

0.9, to be consistent with what was done in the cloud 
parameter retrieval. 


Rejection Criteria 


A number of tests are done to assess the quality of the 
retrieval. The major cause of rejection is difficulty in 
dealing with the effects of clouds on the AIRS radiances. 

1 . Assessment of Cloud Clearing Fit. 


Equations 11 and 10 give the solution for the vector £ 
and the resultant clear column radiances R^. If a 
successful solution is produced, the ensemble R* for 
the cloud clearing channels i should match the 
incoming estimates of clear column radiances R^clr 

to a reasonable degree. A poor match is indicative of 
either a particularly poor first guess or problems in 
handling the effects of clouds on the radiances. We 
compute the weighted residual of the clear column 
radiances used in the computation of £ in brightness 
temperature units 


vl/2 


AF: 


£(V R i,CLR> 2 N ii ! 

J 


INh 1 


/ 3B i 

3T 


\2 


TL 


(59) 


and reject the solution if AF computed when 
generating is greater than 1.75K. Equation 59 is 
equivalent to taking the residual of clear column 
brightness temperatures weighted by the channel noise 
covariance in brightness temperature units. 

2. Difficult Cloud Cases. 


Cases with extensive cloud cover and low contrast are 
particularly difficult. The solution is rejected if the sum 
of the final retrieved cloud fractions for all cloud layers is 
greater than 80% or the total cloud fraction is greater than 
50%, and the total cloud below 500 mb is greater than 
10% and the noise amplification factor is greater than 2, 
or the noise amplification factor is greater than 3, or the 


effective noise amplification factor is greater than 8. We 
also reject cases if the total cloud liquid water determined 
by the microwave product is greater than 0.03 gm/cm\ 

3. Large Residuals in Second Pass Retrievals. 

The general iterative solution is terminated when either 
the residual R n (Equation 44) is less than 10% of_the 
RSS of the predicted noise for each mode AB,, 

(Equation 38) or R° is more than 75% of R n_1 
Slow convergence may indicate a poor solution. We 
reject the solution ifthe converged value of R is greater 
than die RSS of 6B ( in either the surface parameter 

retrieval or the temperature profile retrieval in the 
second pass. Poor convergence generally indicates 
problems with the clear column radiances R , 4 . 

4. Inconsistency of Test “Microwave Only” 
and Combined IR/Microwave Retrievals. 

Under some conditions, the clear column radiances R; 4 
can be poor but all convergence tests are passed. 
Nevertheless, the test microwave only retrieval will 
produce low level temperatures which differ 
significantly from those of the second pass retrieval. 
This generally indicates poor clear column radiances. 
The solution is rejected if tire RMS differences between 
the temperature in the lowest 3 km of the test 
microwave only retrieval differs from that of the second 
pass retrieval by more than 1 .25K. 


4. SIMULATION STUDY 

The simulation study is based on radiances computed 
from conditions derived from a global simulation using 
a version of the operational general circulation model 
(GCM) from NOAA NCEP for December 15, 2000 
(Juang, 1997). Details of the methodology to simulate 
the surface and atmospheric conditions for each AIRS 
footprint are given in Fishbein et al„ 2002. 

The test set is the first scan line of every granule (6 
minute period) for December 15, 2000. The dependent 
data set, on which the regression coefficients are based, 
is taken as cloud free radiances computed from the 
whole day December 10, 2000. A first guess of all the 
geophysical parameters (including surface pressure) 
taken from an 18 hour GCM forecast is available for 
use in the retrieval. We used only the forecast surface 
pressure as it is felt that use of a model forecast first 
guess temperature-moisture profile is unnecessary to 
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analyze AIRS/AMSU data because of the high 
information content of the radiances. 


be assimilated as compared to only 100% clear 
situations. 


The AIRS orbital dataset has the following salient 
features within a given scene made up of nine FOV’s: 

• variable surface topography and surface 
pressure, P, 

• daytime and nighttime conditions 

• temperature, T(P), moisture, q(P), ozone, 
0 3 (P), and other trace constituents from the 
surface to 0.005 mb 

• cloud liquid water profiles, f(P),(only 
affects microwave) 

• multiple level cloud conditions within a 
FOV, with spectrally varying cloud 
emissivity, £ c](1 (u) , and reflectivity, 
pcij(o) , consistent with atmospheric 
conditions. The cloud top pressure, 
emissivity, and reflectivity are spatially 
varying as well. 

• variable surface skin temperature, T,, 
spectral surface emissivity, e(o)and spectral 
surface bi-directional reflectance, p(o) . 

• variable land fraction, with coastlines, 
islands, lakes, etc. 

• orbital simulation with simulated scan lines 
with variable viewing angle and solar 
zenith angle. 

5. RESULTS 

There were a total of 7200 cases in the simulation. In 
74 cases, the microwave retrieval step failed, and no 
retrieval was attempted. Of the remainder, 4604 cases 
were accepted. Figure 2 shows the number of cases, 
and percent accepted, as a function of fractional cloud 
cover, in 0.5% bins. We also show statistics for cases 
we classify as “essentially clear” based on the observed 
radiances. This “essentially clear” flag can be of use to 
the data assimilation community, in which it is common 
to assimilate observed channel radiances under clear 
conditions. The data assimilation community avoids 
cloud contaminated radiances out of fear that noise due 
to handling effects of clouds on the radiances may 
degrade the resultant analyzed Fields. Two potential 
problems with this approach are that only a small 
number of cases are completely clear, thus limiting the 
utility of the sounding data in improvement of forecast 
skill, and scenes with small amounts of cloud cover 
may be mistakenly classified as clear. Our “essentially 
clear” flag is designed to include cases of very small 
amounts of cloudiness to increase the yield of cases to 


A case is a candidate to be called essentially clear if the 
largest eigenvalue of AR'N _1 AR is less than 125 for 
ocean cases and 225 for land cases. This indicates a 
small amount of variability in the radiances in the 3x3 
array of AIRS spots. A larger value of radiance 
variability is allowed for land cases to be called clear 
because surface variability is larger over land than 
ocean. We also define a cloud correction value, ABT, 
as the average difference of the reconstructed clear 
column brightness temperature, 0i , and the 9 spot 
average brightness temperature, @i , for all channels in 
the window regions between 800 cm 1 and 900 cm 1 . 
For the scene to be declared “essentially clear”, ABT 
must be less than or equal to 0.1K. In addition, the 
retrieval must be accepted. 431 cases were called 
essentially clear. Figure 2 includes the number of cases 
called “essentially clear” as a function of actual cloud 
cover. The average cloudiness of all cases was 37.14%, 
t he average cloudiness of all accepted cases was 
31.31%, and the average cloudiness of all cases called 
essentially clear was 0.94%. The percentage of cases 
accepted drops slowly with increasing cloud cover 


Percent Yield vs. Cloud Fraction 



Figure 2. Number of cases as a function of cloud 
fraction (black), percentage of successful cases 
(mauve), and percentage of cases called essentially 
clear (blue). 

until about 50%, and then drops off more rapidly after 
that. Roughly 40% of those accepted cases with cloud 
cover less than 0.5% were identified as “essentially 
clear.” 

Figure 3 shows statistics for the clear column brightness 
temperatures for the 431 essentially clear cases, with 
biases shown in Figure 3a and RMS values shown in 
Figure 3b. The top panel shows the difference between 
the noise free brightness temperatures computed from 
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the truth for a given scene and the average of the 
observed brightness temperatures in the 3x3 array of 
AIRS spots in the scene. This is the correction needed 
to make the observed brightness temperatures match the 
true values. The second panel shows the difference 
between the reconstructed brightness temperatures and 
the average observed values. This is the correction 
made in the cloud clearing. The third panel shows the 
difference between the reconstructed clear column 
brightness temperature and that computed from the 
truth. This is the error in the reconstructed clear 
column brightness temperature. Also shown in the third 
panel of 

the RMS statistics is the single spot channel noise. 

In the mean sense, “essentially clear” spots needed an 
average cloud correction of roughly 0.1K in the 800 
cm 1 - 1150 cm 1 region, and essentially none was made 
on the average. This resulted in a small cold bias in 
this window region in the reconstructed clear column 
brightness temperatures. In the RMS sense, corrections 
of up to 0.25K were needed in the long wave window 
for these cases (some of this is due to channel noise) 
and corrections of about 0.1K were made. For the 


Clear Co I jmn Brightness Temper ature Error 
Essentially Clear Cases (1-31/7200 ) 

0 I AS 



Figure 3a. Mean value of cloud correction needed, 
cloud correction made, and errors of cloud cleared 
brightness temperature for essentially clear cases. 


Clear Column Brightness Tempe r a t ur e Error 
Essentially Clear Cases (431/7200) 
RMS 



Figure 3b. RMS values of cloud correction needed, 
cloud correction made, and cloud cleared brightness 
temperature errors for essentially clear spots. Single 
spot noise is also shown. 

most part, the RMS values of the reconstructed 
brightness temperatures were comparable to, or smaller 
than, the single spot channel noise. Lower values can 
arise if either the channel is considered not to see 
clouds (the noise amplification factor is 1/3) or the 
scene is considered clear or contains very small values 
of t| , resulting in noise amplification factors less than 
1, provided accurate values of t| are obtained. 
Radiances for “essentially clear” cases are definitely 
suitable for data assimilation purposes. 

Figure 4 shows analogous statistics for the 4604 
accepted cases for all cloud conditions. On the average, 
cloud corrections of almost 12K were needed in the 
longwave window region, and the correction made was 
slightly smaller than needed, with about a 0.5K 
negative bias in reconstructed clear column brightness 
temperatures at the worst frequencies. 


Clear Column Brightness Temperature Error 
All Accepted Cases (4604/7200) 

0 t AS 



Figure 4a. Mean values of cloud correction needed, 
cloud correction made, and errors of cloud cleared 
brightness temperatures for all accepted cases. 


Clear Column Brightness Temperature Error 
All Accepted Cases (4604/7200) 

RMS 
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Figure 4b. RMS values of cloud correction needed, 
cloud correction made, and cloud cleared brightness 
temperature errors for all cases. Single spot noise is 
also shown. 

In the RMS sense, reconstructed clear column 
brightness temperatures were still comparable to 
channel noise throughout most of the temperature 
profile sounding regions (650 cm 1 - 750 cm 1 and 
2200 cm 1 - 2400 cm 1 ), but larger than the noise 
elsewhere in the spectrum. RMS errors in the water 
vapor sounding region are still very small and 
radiances in these channels, as well as those in the 
temperature sounding region, should be suitable for 
data assimilation. We encourage researchers in the 
field of data assimilation to test the use of radiances 
for all accepted cases. This would substantially 
increase the number of cases which can be used and 
should further improve forecast skill compared to use 
of radiances in just clear or essentially clear cases. 

Figure 5 shows RMS temperature errors for the 4604 
accepted cases, as well as for the 431 essentially 
clear cases. * Errors are shown for layer mean 
temperatures in roughly 1 km layers from the surface 
to 300 mb, 3 km layers from 300 mb to 30 mb, and 5 
km layers from 30 mb to 1 mb. Results are shown for 
both the regression guess and the final physical 
retrieval. Errors of the surface skin temperature are 
indicated in the figure, as well as average RMS 
temperature profile errors over the layers 700 mb to 
the surface and 100 mb to the surface. The physical 
retrieval improves considerably over the regression 
based retrieval in both clear and cloudy cases with 
the largest improvement near the surface. In cloudy 
cases, part of this improvement is due to use of more 


LAYER MEAN RMS TEMPERATURE ERRORS (°C) 
All Accepted Cases (53.9%) 



TEMP trop 700-surf 

1.45 0.97 1.31 All Coses Regression 

0.59 a. 82 0.92 All Coses Retrieval 

0.50 0.77 0.88 Clear Cases - Regression 

0.25 0.67 0.65 Clear Cases Retrieval 

Figure 5. RMS temperature profile errors. 

accurate clear column radiances as the final physical 
retrieval uses r| 4 while the regression uses p 1 . This 
is not a factor in clear cases however, which still 
show significant improvement in RMS errors near the 
surface. Retrievals under the multi layer cloud cover 
used in this simulation (average cloudiness of 
the accepted cases is 31%) degrade over those in 
clear situations beneath 150 mb, but are still of high 
accuracy. Average tropospheric RMS errors are 
0.82K, and average lower tropospheric errors are 
0.92K, both exceeding or essentially meeting the IK 
RMS error requirement for AIRS. The RMS error in 
the lowest 1 km (1.06K) slightly exceeds IK 
however. Skin temperature errors are 0.25K for clear 
cases and 0.59K for all cases. These include land 
cases and are affected by uncertainties in surface 
spectral emissivity. 

Figure 6 shows mean and RMS errors of retrieved 
surface skin temperature and 1 km tropospheric 
layers up to 344 mb for accepted retrievals as a 
function of actual cloud cover. Also shown is the 
mean and RMS clear column brightness temperature 
error for the 937.8 cm 1 window channel. Mean and 
RMS errors for most temperatures are not very 
sensitive to cloud fraction. Negative biases are 
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found, which increase slowly with increasing cloud 
cover, for the long wave window channel radiance, 
the surface skin temperature, and the temperature in 
the lowest km of the atmosphere. RMS errors of all 
parameters increase slowly with increasing cloud 
cover, especially for the lowest km of the 
atmosphere. 


AIRS Me on Temp e r n tu r e Errors vs. Cloud Fraction 
December 15. 2000 Qata 
All Accepted Cases 



Figure 6a. Mean temperature errors as a function of 
cloud cover for accepted retrievals. 


AIRS RMS T emp a r □ t u r e Errors vs. Cloud Fraction 
December 16. 2000 Data 

Al I Accepted Coses 



Figure 6b. RMS temperature errors as a function of 
cloud cover for accepted retrievals. 

Figure 7 shows RMS percent errors of the retrievals, 
weighted by water vapor amount in the layer, for 
integrated column water vapor in roughly 1 km layers 
from the surface to 200 mb. The RMS percent errors 
of total precipitable water are also indicated in the 
figure. Clouds do not degrade the retrieval profile 
accuracy appreciably. Part of this result may be due 
to sampling differences between clear and cloudy 
areas. 


1 Km LAYER PRECIPITABLE WATER PERCENT ERRORS 
All Accepted Coses (63.9%) 
Essentially Clear Cases (6.0%) 



9.0 16.7 AM Cases Regression 

3,8 16.7 All Cases Retrieval 

10.4 20,6 Clear Coses Regression 

1.3 18.5 Clear Coses — Retrieval 

Figure 7. RMS humidity profile % errors. 

Clouds degrade the total precipitable water accuracy 
by 2.5%, but the 3.8% error of total precipitable 
water for all cases is still extremely good. 1 km layer 
precipitable water errors in the troposphere are 
generally better than 20% in the cloudy cases in 
layers up to about 235 mb. Clear cases RMS % errors 
appear worse than all case errors in the upper 
troposphere. This is probably the result of clear cases 
being considerably dryer than cloudy cases and more 
difficult to retrieve on a percentage basis. 

Figure 8 shows weighted percent errors in the ozone 
profile retrievals in roughly 4 km layers from 260 mb 
to 2.15 mb and in one coarse layer from 260 mb to 
the surface. Also shown is the percentage error in 
total ozone, which is 2.4% for clear cases and 2.6% 
for cloudy cases. The RMS profile errors are better 
than 8% in all layers in both clear and cloudy cases. 
The physical retrieval improves tropospheric ozone 
retrievals considerably over what is obtained by 
regression. 
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Ozone Profile RMS % Errors 
AM Accepted Cases (6.3.9%) 



4.C All Cas«3 Regression 

2.6 All Cases Retrieval 

3.3 Clear Cases Regression 

2.4 Clear Cases - Retrieval 

Figure 8. Ozone profile retrieval errors. 

The accuracy of cloud parameters cannot be compared 
in a straightforward manner because of the existence 
of two cloud layers in some (most) scenes. The 
retrieved cloud fractions are effective, both because of 
errors in retrieved cloud top pressure and assumed 
cloud spectral emissivity. In the retrieval process, 
cloud spectral emissivity was always assumed to be 
0.9 while the true cloud spectral emissivity varied a 
few percent from that. Two straightforward 
parameters to compare are total cloud fraction ft] + 
a 2 > and OLR. The OLR validates 

the cloud parameters, as well as all other parameters, 
in a radiative sense. Clear sky OLR (Fclr) is also 
useful to validate all parameters with the exception of 
clouds. Accuracies should be better for accepted 
cases than for rejected cases, because cloud products 
and all other parameters for rejected cases are based 
on the AMSU retrievals which are less accurate. 

Table 4 shows statistics for retrieved cloud cover, 
OLR and clear sky OLR for accepted and rejected 
cases. The errors for rejected retrievals are poorer 
than for accepted retrievals, but all products should 
be useful for climate studies. The OLR product is 
complementary to OLR measured more directly from 
CERES, also on the Aqua platform, in that the AIRS 


derived product will explain variations of OLR and 
clear sky OLR in space and time in terms of 
variations of surface and atmospheric parameters, 
including cloud cover and height. 


Table 4. Retrieved Cloud Fraction and OLR Errors 


Accepted Cases 

Rejected Cases 

Number 

4604 

2522 

Average Cloud Cover 

31.31% 

47.12% 

Bias 

1.98% 

-1.17% 

RMS Error 

6.33% 

11.75% 

Average OLR 

221.8 W/m 2 

196.2 W/m 2 

Bias 

-1.28 W/m 2 

0.38 W/m 2 

RMS Error 

2.94 W/m 2 

5.20 W/m 2 

Average CLR OLR 

253.0 W/m 2 

238.8 W/m 2 

Bias 

-1.60 W/m 2 

-1.81 W/m 2 

RMS Error 

2.57 W/m 2 

6.76 W/m 2 
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